In [4]:
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.backends.backend_pdf import PdfPages
In [5]:
%matplotlib inline
I am adding white noise to the sine wave to simulate a more realistic wave. Then I am plotting the sine waves with different wave numbers , the filtered wave and the difference between filtered and original wave.
In [6]:
def medianFilter( data, windowLength ):
if (windowLength < len(data)and data.ndim == 1):
tempret = np.zeros(len(data)-windowLength+1) # creating an array where the filtered values will be saved in
if windowLength % 2 ==0: # check if the window length is odd or even because with even window length we get an unsynchrone
for c in range(0, len(tempret)):
tempret[c] = np.median( data[ c : c + windowLength +1 ] ) # write the values of the median filtered wave in tempret, calculate the median of all values in the window
return tempret
else:
for c in range(0, len(tempret)):
tempret[c] = np.median( data[ c : c + windowLength ] )
return tempret
else:
raise ValueError("windowLength must be smaller than len(data) and data must be a 1D array")
In [7]:
def medianSinPlot( waveNumber, windowLength ):
data = np.fromfunction( lambda x: np.sin((x-windowLength / 2)/128 * 2 * np.pi * waveNumber), (128 + windowLength / 2, ) ) #creating an array with a sine wave
noise = np.random.normal(0,0.2,(128 + windowLength / 2)) # creating the noise as an array, filled with random numbers, with the same length as the data array
signal = data + noise # generate the noised signal
datafiltered = medianFilter(signal, windowLength) #calculate the filtered wave with the medianFiltered function
signal = signal[ windowLength / 2 : - windowLength ] # slice the data array to synchronize both waves
datafiltered = datafiltered[ : len(signal) ] # cut the filtered wave to the same length as the data wave
plt.plot( signal )
plt.plot( datafiltered )
plt.plot( signal-datafiltered )
In [8]:
pp = PdfPages( 'median sin with white noise.pdf')
for z in range (2,8):
fig = plt.figure(z, figsize=(30,20)) #creating different figures in one plot, z is the window length
for x in range(1, 5):
for y in range(1, 6):
plt.subplot(5, 5, x + (y-1)*4) #creating different subplots in one figure, with x and y the wave number is calculated
wavenum = (x-1) + (y-1)*4
medianSinPlot( wavenum, z )
plt.suptitle('Median filtered, noised sine waves with window length ' + str(z), fontsize = 60)
plt.xlabel(("Wave number = "+str((x-1) + (y-1)*4)), fontsize=18)
pp.savefig(fig)
pp.close()
In [ ]: